Paso 5 (0)- Trayectorias de hospitalización y mortalidad con foco en condiciones vinculadas a trastornos de salud mental y consumo de sustancias posterior a un primer ingreso por alguno de estos trastornos, en usuarios/as jóvenes y adultos emergentes de población general y pertenecientes a pueblos originarios, 2018-2021, Chile
Preparar la base de datos para representar trayectorias de hospitalización. Mostrar las pruebas de permutaciones (R=1000) para la estabilidad de las soluciones de conglomerados en términos de los índices de calidad.
Autor/a
Andrés González Santa Cruz
Fecha de publicación
1 de abr, 2025
Configurar
Código
# remover objetos y memoria utilizadarm(list=ls());gc()
used (Mb) gc trigger (Mb) max used (Mb)
Ncells 598819 32.0 1294971 69.2 686845 36.7
Vcells 1145005 8.8 8388608 64.0 1879001 14.4
#elegir repositorioif(Sys.info()["sysname"]=="Windows"){options(repos =c(CRAN ="https://cran.dcc.uchile.cl/"))}options(install.packages.check.source ="yes") # Chequea la fuente de los paquetes#borrar caché#system("fc-cache -f -v")if(!require(pacman)){install.packages("pacman");require(pacman)}pacman::p_unlock(lib.loc =.libPaths()) #para no tener problemas reinstalando paquetesif(Sys.info()["sysname"]=="Windows"){if (getRversion() !="4.4.1") { stop("Requiere versión de R 4.4.1. Actual: ", getRversion()) }}cat("quarto version: "); system("quarto --version")
seq_mean_t_dos_grupos <-function(bd =NULL, group1, group2) {# Agrupar por ambas variables resultados <-by(bd, list(group1, group2), seqmeant)# Obtener todas las combinaciones posibles de los grupos combinaciones <-expand.grid(group1 =unique(group1), group2 =unique(group2), stringsAsFactors =FALSE)# Extraer los resultados y asociarlos con las combinaciones resultados_df <-do.call(rbind, lapply(seq_along(resultados), function(i) { group_name1 <-attr(resultados, "dimnames")[[1]][i] group_name2 <-attr(resultados, "dimnames")[[2]][i]data.frame(factor_inclusivo_1 = group_name1, factor_inclusivo_2 = group_name2, Mean = resultados[[i]]) }))# Unir los resultados con las combinaciones para rellenar los valores faltantes final_df <-merge(combinaciones, resultados_df, by.x =c("group1", "group2"), by.y =c("factor_inclusivo_1", "factor_inclusivo_2"), all.x =TRUE)return(final_df)}multinom_pivot_wider <-function(x) {# check inputs match expectatations# create tibble of results df <- tibble::tibble(outcome_level =unique(x$table_body$groupname_col)) df$tbl <- purrr::map( df$outcome_level,function(lvl) { gtsummary::modify_table_body( x, ~dplyr::filter(.x, .data$groupname_col %in% lvl) |> dplyr::ungroup() |> dplyr::select(-.data$groupname_col) ) } )tbl_merge(df$tbl, tab_spanner =paste0("**", df$outcome_level, "**"))}best_subset_multinom <-function(y, x.vars, data) {# y Nombre de la variable dependiente (cadena de texto)# x.vars Vector de nombres de predictores (caracter)# data Dataframe con los datos de entrenamiento# Cargar las librerías necesariasrequire(dplyr)require(purrr)require(tidyr)require(nnet)require(MASS)# Generar todas las combinaciones posibles de predictores predictors_list <-lapply(1:length(x.vars), function(i) {combn(x.vars, i, simplify =FALSE) }) |>unlist(recursive =FALSE)# Inicializar una lista para almacenar los resultados results <-list()# Iterar sobre cada combinación de predictoresfor (i inseq_along(predictors_list)) { predictors <- predictors_list[[i]] formula <-as.formula(paste(y, "~", paste(predictors, collapse ="+")))# Ajustar el modelo multinomial model <-tryCatch( nnet::multinom(formula, data = data, trace =FALSE),error =function(e) NULL )# Si el modelo se ajustó correctamente, almacenar los resultadosif (!is.null(model)) {# Extraer el AIC del modelo aic <-AIC(model)# Almacenar la información en una lista results[[length(results) +1]] <-list(predictors = predictors,model = model,AIC = aic ) } }# Convertir la lista de resultados en un dataframe results_df <- results |> purrr::map_df(function(res) {data.frame(predictors =paste(res$predictors, collapse ="+"),AIC = res$AIC,stringsAsFactors =FALSE ) })# Ordenar los modelos por AIC de menor a mayor results_df <- results_df |>arrange(AIC)return(results_df)}best_subset_multinom_interactions <-function(y, x.vars, data) {# y Nombre de la variable dependiente (cadena de texto)# x.vars Vector de nombres de predictores (caracter)# data Dataframe con los datos de entrenamiento# Cargar las librerías necesariasrequire(dplyr)require(purrr)require(tidyr)require(nnet)require(MASS)# Generar todas las combinaciones posibles de predictores (efectos principales) main_effects_list <-lapply(1:length(x.vars), function(i) {combn(x.vars, i, simplify =FALSE) }) |>unlist(recursive =FALSE)# Inicializar una lista para almacenar los resultados results <-list()# Iterar sobre cada combinación de efectos principalesfor (main_effects in main_effects_list) {# Generar términos de interacción de hasta 3 variables interaction_terms <-list()# Para interacciones de 2 variablesif (length(main_effects) >=2) { interaction_terms_2way <-combn(main_effects, 2, function(x) paste(x, collapse =":")) interaction_terms <-c(interaction_terms, interaction_terms_2way) }# Para interacciones de 3 variablesif (length(main_effects) >=3) { interaction_terms_3way <-combn(main_effects, 3, function(x) paste(x, collapse =":")) interaction_terms <-c(interaction_terms, interaction_terms_3way) }# Combinar efectos principales e interacciones all_terms <-c(main_effects, interaction_terms)# Generar todas las combinaciones posibles de términos (incluyendo interacciones)# Solo se incluyen interacciones si sus efectos principales están presentes term_combinations <-list()# Obtener todos los subconjuntos de efectos principales main_effects_subsets <-lapply(1:length(main_effects), function(i) {combn(main_effects, i, simplify =FALSE) }) |>unlist(recursive =FALSE)# Para cada subconjunto de efectos principalesfor (me in main_effects_subsets) {# Iniciar con los efectos principales terms <- me# Incluir interacciones solo si todos sus efectos principales están incluidos possible_interactions <- interaction_terms[sapply(interaction_terms, function(x) { vars_in_interaction <-unlist(strsplit(x, ":"))all(vars_in_interaction %in% me) }) ]# Generar todas las combinaciones de interacciones para incluir interaction_subsets <-list(NULL)if (length(possible_interactions) >0) { interaction_subsets <-lapply(1:length(possible_interactions), function(i) {combn(possible_interactions, i, simplify =FALSE) }) |>unlist(recursive =FALSE) }# Para cada combinación de interacciones, crear el conjunto completo de términosfor (ints in interaction_subsets) {if (is.null(ints)) { full_terms <- terms } else { full_terms <-c(terms, ints) }# Añadir a la lista de combinaciones de términos term_combinations <-append(term_combinations, list(full_terms)) } }# Ajustar modelos para cada combinación de términosfor (terms in term_combinations) { formula <-as.formula(paste(y, "~", paste(terms, collapse ="+")))# Ajustar el modelo multinomial model <-tryCatch( nnet::multinom(formula, data = data, trace =FALSE),error =function(e) NULL,warning =function(w) NULL )# Si el modelo se ajustó correctamente, almacenar los resultadosif (!is.null(model)) {# Extraer el BIC del modelo bic <-BIC(model)# Almacenar la información en la lista de resultados results[[length(results) +1]] <-list(predictors =paste(terms, collapse =" + "),model = model,BIC = bic ) } } }# Convertir la lista de resultados en un dataframe results_df <- results |> purrr::map_df(function(res) {data.frame(predictors = res$predictors,BIC = res$BIC,stringsAsFactors =FALSE ) })# Ordenar los modelos por BIC de menor a mayor results_df <- results_df |>arrange(BIC)return(results_df)}best_subset_multinom_interactions_parallel <-function(y, x.vars, data) {# y Nombre de la variable dependiente (cadena de texto)# x.vars Vector de nombres de predictores (caracter)# data Dataframe con los datos de entrenamiento# Cargar las librerías necesarias dentro de la funciónrequire(dplyr)require(purrr)require(tidyr)require(nnet)require(MASS)require(foreach)require(doParallel)require(progressr)# Iniciar los gestores de progresohandlers(global =TRUE)handlers("txt")# Generar todas las combinaciones posibles de predictores (efectos principales) main_effects_list <-lapply(1:length(x.vars), function(i) {combn(x.vars, i, simplify =FALSE) }) |>unlist(recursive =FALSE)# Inicializar una lista para almacenar las fórmulas de los modelos formulas_list <-list()# Generar todas las fórmulas posibles con interacciones hasta de 3 variablesfor (main_effects in main_effects_list) {# Generar términos de interacción de hasta 3 variables interaction_terms <-character(0) # Aseguramos que es un vector de caracteres# Para interacciones de 2 variablesif (length(main_effects) >=2) { interaction_terms_2way <-combn(main_effects, 2, function(x) paste(x, collapse =":"), simplify =TRUE) interaction_terms <-c(interaction_terms, interaction_terms_2way) }# Para interacciones de 3 variablesif (length(main_effects) >=3) { interaction_terms_3way <-combn(main_effects, 3, function(x) paste(x, collapse =":"), simplify =TRUE) interaction_terms <-c(interaction_terms, interaction_terms_3way) }# Generar todas las combinaciones posibles de efectos principales main_effects_subsets <-lapply(1:length(main_effects), function(i) {combn(main_effects, i, simplify =FALSE) }) |>unlist(recursive =FALSE)# Para cada subconjunto de efectos principalesfor (me in main_effects_subsets) {# Iniciar con los efectos principales terms <- me# Identificar interacciones cuyos efectos principales están en 'me'if (length(interaction_terms) >0) { possible_interactions <- interaction_terms[vapply(interaction_terms, function(x) { vars_in_interaction <-unlist(strsplit(x, ":"))all(vars_in_interaction %in% me) }, FUN.VALUE =logical(1)) ] } else { possible_interactions <-character(0) }# Generar todas las combinaciones posibles de estas interacciones interaction_subsets <-list(character(0)) # Incluir el caso sin interaccionesif (length(possible_interactions) >0) { interaction_combinations <-lapply(1:length(possible_interactions), function(i) {combn(possible_interactions, i, simplify =FALSE) }) |>unlist(recursive =FALSE) interaction_subsets <-c(interaction_subsets, interaction_combinations) }# Para cada combinación de interaccionesfor (ints in interaction_subsets) { full_terms <-c(terms, ints)# Crear la fórmula del modelo y almacenarla formula_str <-paste(y, "~", paste(full_terms, collapse ="+")) formulas_list <-append(formulas_list, list(formula_str)) } } }# Eliminar posibles duplicados de fórmulas formulas_list <-unique(formulas_list)# Total de modelos a ajustar total_models <-length(formulas_list)# Iniciar el progreso p <-progressor(steps = total_models)# Ajustar los modelos en paralelo usando foreach results_list <-foreach(i =1:total_models, .packages =c("nnet", "MASS"), .combine ='rbind') %dopar% { formula_str <- formulas_list[[i]] formula <-as.formula(formula_str)# Ajustar el modelo model <-tryCatch( nnet::multinom(formula, data = data, trace =FALSE),error =function(e) NULL,warning =function(w) NULL )# Actualizar el progresop(sprintf("Ajustando modelo %d de %d", i, total_models))# Si el modelo se ajustó correctamente, almacenar los resultadosif (!is.null(model)) { bic <-BIC(model)data.frame(predictors = formula_str,BIC = bic,stringsAsFactors =FALSE ) } else {NULL } }# Convertir los resultados a dataframe y ordenar por BIC results_df <-as.data.frame(results_list) results_df <- results_df |>arrange(BIC)return(results_df)}num_cores <- parallel::detectCores() -1cl <-makeCluster(num_cores)registerDoParallel(cl)#pacman job kableExtra tidyverse cluster WeightedCluster devtools TraMineR TraMineRextras NbClust haven ggseqplot gridExtra Tmisc factoextra reticulate withr rmarkdown quartooptions(knitr.kable.NA ='')#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:##:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#knitr::knit_hooks$set(time_it =local({ now <-NULLfunction(before, options) {if (before) {# record the current time before each chunk now <<-Sys.time() } else {# calculate the time difference after a chunk res <-ifelse(difftime(Sys.time(), now)>(60^2),difftime(Sys.time(), now)/(60^2),difftime(Sys.time(), now)/(60^1))# return a character string to show the time x<-ifelse(difftime(Sys.time(), now)>(60^2),paste("Tiempo que demora esta sección:", round(res,1), "horas"),paste("Tiempo que demora esta sección:", round(res,1), "minutos"))paste('<div class="message">', gsub('##', '\n', x),'</div>', sep ='\n') } }}))knitr::opts_chunk$set(time_it =TRUE)#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:format_cells <-function(df, rows ,cols, value =c("italics", "bold", "strikethrough")){# select the correct markup# one * for italics, two ** for bold map <-setNames(c("*", "**", "~~"), c("italics", "bold", "strikethrough")) markup <- map[value] for (r in rows){for(c in cols){# Make sure values are not factors df[[c]] <-as.character( df[[c]])# Update formatting df[r, c] <-ifelse(nchar(df[r, c])==0,"",paste0(markup, gsub(" ", "", df[r, c]), markup)) } }return(df)}#To produce line breaks in messages and warningsknitr::knit_hooks$set(error =function(x, options) {paste('\n\n<div class="alert alert-danger" style="font-size: small !important;">',gsub('##', '\n', gsub('^##\ Error', '**Error**', x)),'</div>', sep ='\n') },warning =function(x, options) {paste('\n\n<div class="alert alert-warning" style="font-size: small !important;">',gsub('##', '\n', gsub('^##\ Warning:', '**Warning**', x)),'</div>', sep ='\n') },message =function(x, options) {paste('<div class="message" style="font-size: small !important;">',gsub('##', '\n', x),'</div>', sep ='\n') })#_#_#_#_#_#_#_#_#_#_#_#_#_invisible("Function to format CreateTableOne into a database")as.data.frame.TableOne <-function(x, ...) {capture.output(print(x,showAllLevels =TRUE, varLabels = T,...) -> x) y <-as.data.frame(x) y$characteristic <- dplyr::na_if(rownames(x), "") y <- y |>fill(characteristic, .direction ="down") |> dplyr::select(characteristic, everything())rownames(y) <-NULL y}#_#_#_#_#_#_#_#_#_#_#_#_#_# Austin, P. C. (2009). The Relative Ability of Different Propensity # Score Methods to Balance Measured Covariates Between # Treated and Untreated Subjects in Observational Studies. Medical # Decision Making. https://doi.org/10.1177/0272989X09341755smd_bin <-function(x,y){ z <- x*(1-x) t <- y*(1-y) k <-sum(z,t) l <- k/2return((x-y)/sqrt(l))}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:format_table_vec <-function(tbl, digits =1) { counts <-as.numeric(tbl) percentages <-prop.table(tbl) *100 formatted <-sapply(seq_along(counts), function(i) { p_val <- percentages[i]# Si el porcentaje es prácticamente entero, formatea sin decimalesif (abs(p_val -round(p_val)) < .Machine$double.eps^0.5) { p_str <-sprintf("%.0f", p_val) } else { p_str <-sprintf(paste0("%.", digits, "f"), p_val) }paste0(counts[i], " (", p_str, ")") }) formatted}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:if(.Platform$OS.type =="windows") withAutoprint({memory.size()memory.size(TRUE)memory.limit()})
> memory.size()
[1] Inf
> memory.size(TRUE)
[1] Inf
> memory.limit()
[1] Inf
Código
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:func_tab_range_clus<-function(range_clus){rbind.data.frame(lapply(list(as.vector(rev(sort(table(range_clus$clustering$cluster2)))),as.vector(rev(sort(table(range_clus$clustering$cluster3)))),as.vector(rev(sort(table(range_clus$clustering$cluster4)))),as.vector(rev(sort(table(range_clus$clustering$cluster5)))),as.vector(rev(sort(table(range_clus$clustering$cluster6)))),as.vector(rev(sort(table(range_clus$clustering$cluster7)))),as.vector(rev(sort(table(range_clus$clustering$cluster8)))),as.vector(rev(sort(table(range_clus$clustering$cluster9)))),as.vector(rev(sort(table(range_clus$clustering$cluster10)))),as.vector(rev(sort(table(range_clus$clustering$cluster11)))),as.vector(rev(sort(table(range_clus$clustering$cluster12)))),as.vector(rev(sort(table(range_clus$clustering$cluster13)))),as.vector(rev(sort(table(range_clus$clustering$cluster14)))),as.vector(rev(sort(table(range_clus$clustering$cluster15)))) ),function(x) { length_out <-max(sapply(list(as.vector(rev(sort(table(range_clus$clustering$cluster2)))),as.vector(rev(sort(table(range_clus$clustering$cluster3)))),as.vector(rev(sort(table(range_clus$clustering$cluster4)))),as.vector(rev(sort(table(range_clus$clustering$cluster5)))),as.vector(rev(sort(table(range_clus$clustering$cluster6)))),as.vector(rev(sort(table(range_clus$clustering$cluster7)))),as.vector(rev(sort(table(range_clus$clustering$cluster8)))),as.vector(rev(sort(table(range_clus$clustering$cluster9)))),as.vector(rev(sort(table(range_clus$clustering$cluster10)))),as.vector(rev(sort(table(range_clus$clustering$cluster11)))),as.vector(rev(sort(table(range_clus$clustering$cluster12)))),as.vector(rev(sort(table(range_clus$clustering$cluster13)))),as.vector(rev(sort(table(range_clus$clustering$cluster14)))),as.vector(rev(sort(table(range_clus$clustering$cluster15)))) ), length))c(x, rep(NA, length_out -length(x))) } ))|>t() |>data.frame()|>`rownames<-`(NULL)}frobenius_norm <-function(matrix1, matrix2) {if (!all(dim(matrix1) ==dim(matrix2))) {stop("Matrices must have the same dimensions") }# Replace NA values with 0 (or any other desired default) matrix1[is.na(matrix1)] <-0 matrix2[is.na(matrix2)] <-0# Calculate the residuals residuals <- matrix1 - matrix2# Frobenius norm frobenius <-sqrt(sum(residuals^2))return(frobenius)}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:print.seqnullcqi.powder <-function(x, norm =FALSE, quant =0.95, digits =2, append =FALSE, ...) { confcqi2 <-function(nullstat, quant, n){ alpha <- (1-quant)/2#calpha <- alpha+(alpha-1)/n#print(c(calpha, alpha))#minmax <- quantile(nullstat, c(calpha, 1-calpha)) minmax <-quantile(nullstat, c(alpha, 1-alpha))return(minmax) } normstatcqi2 <-function(bcq, stat, norm=TRUE){ origstat <- bcq$clustrange$stats[, stat] nullstat <- bcq$stats[[stat]]#normstat <- rbind(nullstat, origstat)if(norm){for(i inseq_along(origstat)){ mx <-mean(nullstat[, i]) sdx <-sd(nullstat[, i]) nullstat[ , i] <- (nullstat[, i]-mx)/sdx origstat[i] <- (origstat[i]-mx)/sdx } } alldatamax <-apply(nullstat, 1, max)#as.vector(xx) sumcqi <-list(origstat=origstat, nullstat=nullstat, alldatamax=alldatamax)return(sumcqi) }cat("Parametric bootstrap cluster analysis validation\n")cat("Sequence analysis null model:", deparse(x$nullmodel), "\n")cat("Number of bootstraps:", x$R, "\n")cat("Clustering method:", ifelse(x$kmedoid, "PAM/K-Medoid", paste0("hclust with ", x$hclust.method)), "\n")cat("Seqdist arguments:", deparse(x$seqdist.args), "\n\n\n") alls <-as.data.frame(x$clustrange$stats) quants <-rep("", ncol(alls))names(quants) <-colnames(alls)for (ss incolnames(alls)) { sumcqi <-normstatcqi2(x, stat = ss, norm = norm) alls[, ss] <-as.character(round(sumcqi$origstat, digits = digits)) borne <-as.character(round(confcqi2(sumcqi$alldatamax, quant, x$R), digits = digits)) quants[ss] <-paste0("[", borne[1], "; ", borne[2], "]") } results_tibble <- tibble::as_tibble(rbind(alls, rep("", length(quants)), quants))# Print a summary to the console for immediate feedbackrownames(results_tibble) <-c(rownames(x$clustrange$stats), "", paste("Null Max-T", quant, "interval")) results_df <-as.data.frame(results_tibble)print(results_tibble, ...)return(list(results_tibble= results_tibble, results_df= results_df ))}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:##:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:# Función para aplicar la prueba de Fisher a todas las combinaciones de filas usando todas las columnasfisher_posthoc_all_cols <-function(contingency_table) {# Obtener combinaciones de filas (pares) row_pairs <-combn(rownames(contingency_table), 2, simplify =FALSE)# Aplicar la prueba de Fisher a cada par de filas usando todas las columnas al mismo tiempo results <-map_dfr(row_pairs, function(pair) {# Crear tabla de 2xN para el par de filas en todas las columnas sub_table <- contingency_table[pair, , drop =FALSE]# Aplicar el test de Fisher test_result <-fisher.test(sub_table, simulate.p.value=T,B=1e4)# Devolver los resultados en un data frametibble(Row1 = pair[1],Row2 = pair[2],p.value = test_result$p.value ) })# Ajustar p-valores usando el método de Holm results <- results |>mutate(p.adjusted =p.adjust(p.value, method ="holm"))return(results)}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:save_base_plot_as_grob <-function(plot_expr, res=300, width =1600, height=1200) {# Crea un archivo temporal con extensión .png filename <-tempfile(fileext =".png")# Guarda el gráfico en alta resolución en el archivo temporalpng(filename, width = width, height = height, res = res)replayPlot(plot_expr) # Reproduce el gráfico grabadodev.off() # Cierra el dispositivo gráfico# Convierte el archivo PNG en un objeto gráfico (grob) grob <- grid::rasterGrob(png::readPNG(filename), interpolate =TRUE)return(grob) # Devuelve el grob}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:chisq_cramerv<-function(contingency_table){ chisq_test <-chisq.test(contingency_table) cramers_v <-sqrt(chisq_test$statistic / (sum(contingency_table) * (min(dim(contingency_table)) -1)))list(chisq_statistic=sprintf("%1.2f", chisq_test$statistic), chisq_df= chisq_test$parameter, chisq_p_value =ifelse(chisq_test$p.value<.001, "<0.001", sprintf("%1.4f", chisq_test$p.value)), cramers_v =sprintf("%1.2f", cramers_v))}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#oneway_anova_effect_size <-function(values, group) {# Perform one-way ANOVA anova_result <-aov(values ~ group)# Summarize ANOVA results anova_summary <-summary(anova_result)# Extract sums of squares ss_between <- anova_summary[[1]]$"Sum Sq"[1] ss_total <-sum(anova_summary[[1]]$"Sum Sq")# Calculate eta-squared eta_squared <- ss_between / ss_total# Return ANOVA summary and effect sizelist(anova_summary = anova_summary,eta_squared = eta_squared )}
head(arrange(data.frame(table(diag_todos)), -Freq) |> dplyr::mutate(perc=scales::percent(Freq/sum(Freq), accuracy=.01)),20) |> knitr::kable("markdown", caption="Diagnósticos más frecuentes, código CIE-10, detalle y frecuencia")
Diagnósticos más frecuentes, código CIE-10, detalle y frecuencia
diag_todos
Freq
perc
F329
1089
4.75%
F322
959
4.19%
F609
846
3.69%
F603
606
2.65%
F432
548
2.39%
F209
449
1.96%
F192
432
1.89%
F319
399
1.74%
Z915
391
1.71%
F200
304
1.33%
F323
258
1.13%
F29X
249
1.09%
F121
222
0.97%
F101
208
0.91%
F419
202
0.88%
F142
185
0.81%
F321
179
0.78%
F102
170
0.74%
F171
167
0.73%
T509
163
0.71%
Código
# **F329 (n=990 a 1089)** – *Episodio depresivo no especificado.* # Se refiere a un cuadro depresivo cuyos síntomas no cumplen criterios completos para especificar la gravedad o características particulares.# **F322 (n=845 a 959)** – *Episodio depresivo grave sin síntomas psicóticos.* # Episodio depresivo intenso que no presenta alucinaciones ni ideas delirantes, pero con afectación significativa del funcionamiento.# **F609 (n=770 a 846)** – *Trastorno de la personalidad sin especificar.* # Diagnóstico que incluye rasgos de personalidad patológicos que no se ajustan a categorías específicas conocidas.# **F603 (n=550 a 606)** – *Trastorno de la personalidad emocionalmente inestable (tipo límite).* # También llamado “trastorno límite de la personalidad”, caracterizado por inestabilidad emocional, relaciones interpersonales conflictivas y conducta impulsiva.# **F432 (n=491 a 548)** – *Trastornos de adaptación.* # Reacciones emocionales y/o conductuales que surgen como respuesta a un cambio o factor estresante identificable, dificultando la adaptación normal.# **F209 (n=433 a 449)** – *Esquizofrenia no especificada.* # Forma de esquizofrenia en la que no se pueden determinar subtipos (paranoide, catatónica, etc.) o faltan detalles para clasificarlos.# **F192 (n=394 a 432)** – *Síndrome de dependencia por uso de múltiples drogas.* # Dependencia y uso problemático de diversas sustancias psicoactivas, con patrones de consumo repetitivo y dificultades para el control.# **F319 (n=369 a 399)** – *Trastorno bipolar no especificado.* # Forma de trastorno bipolar con episodios de alteración del estado de ánimo, donde faltan datos para clasificar un subtipo específico.# **Z915 (n=360 a 391)** – *Antecedentes personales de autolesiones.* # Historia previa de conducta autolesiva o intento de suicidio, utilizada para codificar factores influyentes en el estado de salud actual.# **F200 (n=292 a 304)** – *Esquizofrenia paranoide.* # Subtipo de esquizofrenia caracterizado principalmente por la presencia de delirios y alucinaciones de tipo paranoide.# **F29X (n=232 a 249)** – *Psicosis no orgánica no especificada.* # Trastorno psicótico sin evidencia de causa orgánica, cuyos rasgos no son suficientes para un diagnóstico más preciso.# **F323 (n=230 a 258)** – *Episodio depresivo grave con síntomas psicóticos.* # Episodio depresivo intenso que incluye delirios, alucinaciones u otras manifestaciones psicóticas.# **F121 (n=200 a 222)** – *Uso perjudicial de cannabis.* # Consumo de cannabis que causa un deterioro en el funcionamiento personal o social, sin llegar al síndrome de dependencia.# **F101 (n=187 a 208)** – *Uso perjudicial de alcohol.* # Patrón de consumo de alcohol que provoca daño a la salud física o mental, sin cumplir criterios de dependencia.# **F419 (n=176 a 202)** – *Trastorno de ansiedad no especificado.* # Ansiedad significativa y persistente que no se encuadra en categorías específicas (p.ej. fobias, pánico, etc.).# **F142 (n=164 a 185)** – *Síndrome de dependencia de cocaína.* # Presencia de dependencia a la cocaína, con anhelo intenso y dificultad para controlar o interrumpir el consumo.# **F321 (n=163 a 179)** – *Episodio depresivo moderado.* # Estado depresivo con síntomas clínicamente relevantes, pero de gravedad intermedia entre leve y grave.# **F102 (n=157 a 170)** – *Síndrome de dependencia de alcohol.* # Patrón de dependencia caracterizado por la necesidad imperiosa de beber y dificultad para controlar el consumo.# **F171 (n=167)** – *Uso perjudicial de tabaco.* # Consumo de tabaco que genera consecuencias físicas o mentales negativas, sin cumplir criterios de dependencia.
Tiempo que demora esta sección: 0.1 minutos
El algoritmo PAM busca formar clústeres de manera eficiente, pero puede dar resultados poco generalizables si parte de medoides iniciales poco adecuados. Para evitar esto, generamos una solución de medoides utilizando como punto de partida los clústeres obtenidos previamente con un método jerárquico.
Exportamos .RDS para generar las pruebas de permutaciones.
cqi_quarter<-rbind.data.frame(cbind.data.frame(algo="hac", type="om", time="quarter", k=2:15, corr=F, om_dist_quarter_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="hac", type="lcs", time="quarter", k=2:15, corr=F, lcs_dist_quarter_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="quarter", k=2:15, corr=F, pamRange_quarter_om$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="quarter", k=2:15, corr=T, pamRange_quarter_om2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="quarter", k=2:15, corr=F, pamRange_quarter_lcs$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="quarter", k=2:15, corr=T, pamRange_quarter_lcs2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))))|> dplyr::select(algo, type, time, k, corr, PBC, ASW, HC, HG, R2, R2sq)# round(summary(silhouette(as.integer(om_dist_quarter_c$clustering$cluster2), as.dist(dist_quarter_om)))$clus.avg.widths,2)[attr(rev(sort(table(om_dist_quarter_c$clustering$cluster2))),"names")]#functión para generar el gráficotabs_quarter_clus_sol<-rbind.data.frame(func_tab_range_clus(om_dist_quarter_c),func_tab_range_clus(lcs_dist_quarter_c),func_tab_range_clus(pamRange_quarter_om),func_tab_range_clus(pamRange_quarter_om2),func_tab_range_clus(pamRange_quarter_lcs),func_tab_range_clus(pamRange_quarter_lcs2))#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:##:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:# Inicializamos una lista para almacenar los resultadosresultados_list <-list()# Definimos un rango para los clusters a evaluarcluster_range <-2:15# Definimos los métodos y sus variablesmetodos <-list(hac_om =list(data = om_dist_quarter_c, dist = dist_quarter_om),hac_lcs =list(data = lcs_dist_quarter_c, dist = dist_quarter_lcs),pam_om0 =list(data = pamRange_quarter_om, dist = dist_quarter_om),pam_om1 =list(data = pamRange_quarter_om2, dist = dist_quarter_om),pam_lcs0 =list(data = pamRange_quarter_lcs, dist = dist_quarter_lcs),pam_lcs1 =list(data = pamRange_quarter_lcs2, dist = dist_quarter_lcs))# Número máximo de clusters para definir las columnasmax_clusters <-max(cluster_range)# Iteramos sobre cada métodofor (metodo innames(metodos)) {# Creamos un data frame temporal para cada método metodo_result <-data.frame()# Iteramos sobre cada cluster en el rangofor (cluster in cluster_range) {# Construimos el nombre del cluster dinámicamente cluster_name <-paste0("cluster", cluster)# Intentamos calcular los valores de silhouette silhouette_values <-tryCatch(round(summary(silhouette(as.integer(metodos[[metodo]]$data$clustering[[cluster_name]]), as.dist(metodos[[metodo]]$dist)))$clus.avg.widths[attr(rev(sort(table(metodos[[metodo]]$data$clustering[[cluster_name]]))),"names")], 2),error =function(e) rep(NA, cluster) )# Creamos un vector con las columnas llenando con NA si faltan valores silhouette_full <-c(silhouette_values, rep(NA, max_clusters -length(silhouette_values)))# Creamos un data frame temporal con los resultados para este cluster cluster_result <-data.frame(Metodo = metodo,Cluster = cluster,t(silhouette_full) # Transponemos los valores para que cada uno sea una columna )# Nombramos dinámicamente las columnas de silhouettecolnames(cluster_result)[3:(3+ max_clusters -1)] <-paste0("asw", 1:max_clusters)# Añadimos el resultado del cluster al data frame del método metodo_result <-rbind(metodo_result, cluster_result) }# Agregamos los resultados del método a la lista general resultados_list[[metodo]] <- metodo_result}# Combinamos todos los resultados en un único data frameavs_por_cluster_quarter <-do.call(rbind, resultados_list)# Ordenamos por Método y Clusteravs_por_cluster_quarter <- avs_por_cluster_quarter[order(avs_por_cluster_quarter$Metodo, avs_por_cluster_quarter$Cluster), ]bind_cols(cqi_quarter, tabs_quarter_clus_sol)|> dplyr::mutate(corr= dplyr::case_when(corr==TRUE& algo!="hac"~"1",corr==FALSE& algo!="hac"~"0",T~""), key=paste0(algo,"_",type,corr,"_",k))|>left_join(dplyr::mutate(avs_por_cluster_quarter, key=paste0(Metodo,"_",Cluster)), by="key") |> dplyr::select(-Metodo, -Cluster) |>`rownames<-`(NULL) |> dplyr::mutate(calc=round(PBC*(1/HC)*HG,2)) |> dplyr::arrange(desc(ASW)) |> dplyr::select(c("algo", "type", "time", "k", "corr", "PBC", "ASW", "HC", "HG", "R2", "R2sq", "calc", paste0("X",1:15), paste0("asw",1:15))) |> (\(df) {assign("asw_quarter_qci", dplyr::select(df, -"time"), envir = .GlobalEnv) rio::export(df, "_output/sol_conglomerados_tab_quarter_25.xlsx") knitr::kable(df, "markdown", caption ="CQIs y frecuencias en conglomerados (trimestre)") })()
CQIs y frecuencias en conglomerados (trimestre)
algo
type
time
k
corr
PBC
ASW
HC
HG
R2
R2sq
calc
X1
X2
X3
X4
X5
X6
X7
X8
X9
X10
X11
X12
X13
X14
X15
asw1
asw2
asw3
asw4
asw5
asw6
asw7
asw8
asw9
asw10
asw11
asw12
asw13
asw14
asw15
pam
om
quarter
6
0
0.59
0.59
0.10
0.77
0.38
0.49
4.54
4885
746
358
239
224
174
0.66
0.62
0.22
-0.03
0.49
-0.05
pam
om
quarter
14
0
0.59
0.59
0.06
0.85
0.51
0.59
8.36
4121
745
258
220
169
167
153
147
127
113
107
102
100
97
0.74
0.57
0.23
0.47
0.05
0.26
-0.13
0.06
0.38
-0.06
0.39
0.49
0.34
0.40
pam
om
quarter
15
0
0.60
0.59
0.05
0.86
0.52
0.59
10.32
4042
745
258
220
169
167
153
147
127
111
107
102
100
94
84
0.76
0.57
0.22
0.46
0.04
0.25
-0.13
0.05
0.37
-0.06
0.38
0.48
0.34
0.43
0.32
pam
om
quarter
6
1
0.59
0.59
0.10
0.77
0.38
0.49
4.54
4885
746
358
239
224
174
0.66
0.62
0.22
-0.03
0.49
-0.05
pam
om
quarter
7
1
0.59
0.59
0.10
0.78
0.40
0.51
4.60
4729
746
353
224
212
189
173
0.68
0.61
0.22
0.48
0.04
0.20
-0.05
hac
om
quarter
4
0.50
0.58
0.19
0.66
0.28
0.36
1.74
5352
738
312
224
0.61
0.66
-0.15
0.52
hac
om
quarter
5
0.50
0.58
0.19
0.66
0.30
0.40
1.74
5352
738
277
224
35
0.61
0.66
-0.04
0.52
-0.04
pam
om
quarter
4
0
0.49
0.58
0.20
0.66
0.30
0.37
1.62
5248
750
399
229
0.62
0.62
0.08
0.48
pam
om
quarter
5
0
0.55
0.58
0.14
0.72
0.34
0.45
2.83
5119
746
361
224
176
0.62
0.63
0.23
0.50
-0.04
pam
om
quarter
9
0
0.59
0.58
0.09
0.80
0.44
0.54
5.24
4495
746
351
222
173
169
159
156
155
0.69
0.60
0.18
0.49
0.12
0.33
-0.08
0.21
0.10
pam
om
quarter
10
0
0.59
0.58
0.09
0.81
0.45
0.55
5.31
4380
745
351
221
172
169
159
155
153
121
0.70
0.59
0.17
0.48
0.11
0.32
-0.09
0.09
0.21
0.26
pam
om
quarter
13
0
0.60
0.58
0.06
0.84
0.49
0.58
8.40
4204
745
258
221
169
167
153
147
127
117
116
102
100
0.72
0.58
0.25
0.47
0.06
0.27
-0.13
0.07
0.39
0.30
-0.06
0.51
0.34
pam
om
quarter
4
1
0.49
0.58
0.20
0.66
0.30
0.37
1.62
5248
750
399
229
0.62
0.62
0.08
0.48
pam
om
quarter
5
1
0.55
0.58
0.14
0.72
0.34
0.45
2.83
5119
746
361
224
176
0.62
0.63
0.23
0.50
-0.04
pam
om
quarter
8
1
0.61
0.58
0.08
0.82
0.42
0.52
6.25
4583
746
336
226
212
178
173
172
0.69
0.60
0.21
0.45
0.00
-0.08
0.30
-0.08
pam
om
quarter
13
1
0.60
0.58
0.06
0.84
0.50
0.58
8.40
4204
745
258
221
165
152
144
143
134
128
122
110
100
0.71
0.58
0.25
0.47
0.02
-0.05
0.27
0.09
0.52
0.18
-0.14
0.57
0.34
pam
om
quarter
14
1
0.59
0.58
0.06
0.85
0.51
0.59
8.36
4121
745
258
220
165
152
148
146
128
120
114
107
102
100
0.73
0.57
0.24
0.47
0.00
-0.13
0.05
0.27
0.59
0.12
-0.06
0.60
0.40
0.34
pam
om
quarter
15
1
0.60
0.58
0.05
0.86
0.52
0.59
10.32
4042
745
258
220
165
152
148
146
126
111
106
106
101
100
100
0.75
0.57
0.23
0.46
-0.02
-0.14
0.04
0.22
0.60
-0.05
0.61
0.30
0.11
0.34
0.42
hac
lcs
quarter
15
0.56
0.57
0.02
0.95
0.55
0.65
26.60
3327
894
485
454
257
256
237
203
202
100
85
80
30
10
6
1.00
-0.13
1.00
-0.07
-0.23
0.12
0.16
0.02
0.36
0.02
-0.03
-0.07
0.00
0.19
0.51
pam
om
quarter
7
0
0.59
0.57
0.10
0.78
0.40
0.51
4.60
4729
746
353
225
224
176
173
0.66
0.61
0.22
-0.01
0.48
0.30
-0.06
pam
om
quarter
8
0
0.59
0.57
0.09
0.80
0.42
0.52
5.24
4612
746
351
222
212
169
159
155
0.67
0.60
0.19
0.49
0.01
0.35
-0.07
0.11
pam
om
quarter
11
0
0.58
0.57
0.08
0.81
0.47
0.55
5.87
4308
745
351
221
171
169
159
148
130
121
103
0.69
0.59
0.15
0.48
0.08
0.28
-0.09
0.10
0.37
0.25
0.52
pam
om
quarter
12
0
0.58
0.57
0.08
0.81
0.48
0.56
5.87
4308
745
258
221
171
169
153
147
130
121
103
100
0.69
0.59
0.26
0.48
0.08
0.28
-0.10
0.09
0.37
0.25
0.52
0.34
pam
om
quarter
9
1
0.60
0.57
0.08
0.82
0.44
0.54
6.15
4460
746
336
226
200
174
172
158
154
0.69
0.59
0.19
0.44
-0.03
-0.08
-0.08
0.42
0.22
pam
om
quarter
10
1
0.61
0.57
0.07
0.84
0.46
0.55
7.32
4369
744
350
223
191
160
154
152
142
141
0.71
0.59
0.16
0.46
-0.02
-0.15
0.43
-0.05
0.13
0.32
hac
lcs
quarter
13
0.54
0.56
0.03
0.92
0.53
0.63
16.56
3327
1348
485
257
256
237
203
202
130
85
80
10
6
1.00
-0.16
1.00
-0.23
0.12
0.18
0.03
0.36
-0.08
-0.02
-0.03
0.19
0.54
hac
lcs
quarter
14
0.54
0.56
0.03
0.93
0.54
0.64
16.74
3327
1348
485
257
256
237
203
202
100
85
80
30
10
6
1.00
-0.16
1.00
-0.23
0.12
0.18
0.02
0.36
0.02
-0.03
-0.03
0.00
0.19
0.51
pam
om
quarter
11
1
0.58
0.56
0.08
0.81
0.47
0.55
5.87
4308
745
351
221
168
159
148
134
133
130
129
0.67
0.59
0.16
0.48
0.06
-0.09
0.09
0.26
0.16
0.58
0.40
pam
om
quarter
12
1
0.58
0.56
0.08
0.81
0.48
0.56
5.87
4308
745
258
221
168
153
147
138
130
129
129
100
0.67
0.59
0.27
0.48
0.06
-0.11
0.08
0.24
0.58
0.40
0.17
0.34
hac
om
quarter
3
0.49
0.55
0.20
0.65
0.24
0.31
1.59
5352
962
312
0.62
0.41
-0.15
pam
om
quarter
2
0
0.33
0.55
0.34
0.54
0.18
0.22
0.52
5868
758
0.53
0.63
pam
om
quarter
3
0
0.44
0.55
0.24
0.62
0.24
0.30
1.14
5456
756
414
0.58
0.62
0.07
pam
om
quarter
2
1
0.33
0.55
0.34
0.54
0.18
0.22
0.52
5868
758
0.53
0.63
pam
om
quarter
3
1
0.44
0.55
0.24
0.62
0.24
0.30
1.14
5456
756
414
0.58
0.62
0.07
hac
om
quarter
2
0.36
0.54
0.30
0.55
0.19
0.22
0.66
5664
962
0.56
0.44
pam
lcs
quarter
9
1
0.58
0.54
0.06
0.91
0.42
0.54
8.80
4560
701
302
273
189
186
185
115
115
0.65
0.54
0.31
-0.09
0.16
0.50
0.08
0.03
-0.08
pam
lcs
quarter
11
1
0.58
0.54
0.04
0.93
0.45
0.57
13.48
4385
695
299
224
203
193
190
138
128
95
76
0.67
0.54
0.30
-0.22
0.11
0.04
0.49
0.46
0.16
-0.03
-0.07
hac
lcs
quarter
11
0.54
0.53
0.03
0.92
0.50
0.60
16.56
3327
1348
742
256
237
203
202
130
90
85
6
1.00
-0.16
0.37
0.12
0.18
0.03
0.36
-0.08
-0.09
-0.02
0.54
hac
lcs
quarter
12
0.54
0.53
0.03
0.92
0.51
0.61
16.56
3327
1348
742
256
237
203
202
130
85
80
10
6
1.00
-0.16
0.37
0.12
0.18
0.03
0.36
-0.08
-0.02
-0.03
0.19
0.54
pam
lcs
quarter
7
1
0.56
0.53
0.10
0.87
0.37
0.50
4.87
4758
738
307
239
226
198
160
0.62
0.49
0.29
0.02
0.04
0.47
-0.07
pam
lcs
quarter
8
1
0.57
0.53
0.07
0.89
0.40
0.51
7.25
4604
701
307
274
209
189
187
155
0.64
0.54
0.28
-0.09
0.04
0.17
0.50
-0.07
pam
lcs
quarter
10
1
0.58
0.53
0.05
0.92
0.43
0.55
10.67
4448
701
302
273
189
185
184
145
112
87
0.66
0.53
0.30
-0.11
0.13
0.09
0.52
-0.09
0.43
0.00
pam
lcs
quarter
14
1
0.57
0.53
0.03
0.94
0.49
0.59
17.86
4127
695
299
224
203
193
183
140
128
99
86
85
85
79
0.67
0.51
0.28
-0.29
-0.01
0.02
0.50
0.47
0.12
0.75
0.34
-0.08
0.00
0.72
hac
om
quarter
11
0.55
0.52
0.10
0.79
0.43
0.54
4.35
4210
738
356
344
277
268
224
162
20
15
12
0.69
0.62
0.15
-0.15
-0.16
0.04
0.48
-0.16
0.10
0.07
0.12
hac
om
quarter
12
0.56
0.52
0.09
0.79
0.44
0.55
4.92
4210
738
356
344
268
224
203
162
74
20
15
12
0.69
0.62
0.14
-0.15
0.04
0.48
-0.08
-0.17
-0.05
0.09
0.07
0.12
hac
om
quarter
13
0.56
0.52
0.09
0.79
0.44
0.56
4.92
4210
738
356
344
268
224
203
102
74
60
20
15
12
0.69
0.62
0.14
-0.15
0.04
0.47
-0.08
0.06
-0.06
0.07
0.08
0.01
0.12
hac
lcs
quarter
9
0.51
0.52
0.06
0.86
0.48
0.57
7.31
3327
1807
742
237
202
130
90
85
6
1.00
-0.17
0.38
0.21
0.37
-0.04
-0.03
0.01
0.54
hac
lcs
quarter
10
0.53
0.52
0.04
0.89
0.49
0.58
11.79
3327
1604
742
237
203
202
130
90
85
6
1.00
-0.16
0.38
0.19
0.03
0.37
-0.07
-0.04
-0.02
0.54
pam
lcs
quarter
12
0
0.57
0.52
0.04
0.93
0.46
0.57
13.25
4264
742
294
251
201
191
146
126
125
117
85
84
0.65
0.42
-0.08
-0.18
0.62
0.47
-0.09
0.43
0.67
0.36
0.38
0.01
pam
lcs
quarter
13
0
0.57
0.52
0.04
0.93
0.48
0.57
13.25
4185
742
294
251
201
191
146
133
123
100
92
84
84
0.66
0.40
-0.08
-0.21
0.62
0.46
-0.09
0.55
0.25
0.74
0.53
0.37
0.01
pam
lcs
quarter
6
1
0.53
0.52
0.13
0.83
0.35
0.45
3.38
4821
739
349
281
241
195
0.61
0.50
0.21
-0.06
0.01
0.49
pam
lcs
quarter
12
1
0.57
0.52
0.04
0.92
0.46
0.57
13.11
4286
692
300
299
188
183
138
132
125
109
98
76
0.64
0.53
-0.20
0.29
0.49
0.05
0.28
0.14
0.67
0.44
-0.05
-0.07
pam
lcs
quarter
13
1
0.58
0.52
0.03
0.94
0.48
0.59
18.17
4203
706
299
219
218
192
182
131
110
100
98
92
76
0.67
0.48
0.29
-0.05
-0.13
0.03
0.51
0.12
0.65
0.53
-0.05
0.28
-0.07
hac
om
quarter
7
0.52
0.51
0.16
0.69
0.35
0.47
2.24
4984
738
356
277
224
35
12
0.55
0.65
0.23
-0.11
0.51
-0.06
0.18
hac
om
quarter
8
0.53
0.51
0.13
0.71
0.38
0.49
2.89
4640
738
356
344
277
224
35
12
0.60
0.64
0.20
-0.09
-0.12
0.50
-0.06
0.12
hac
om
quarter
9
0.53
0.51
0.13
0.73
0.40
0.51
2.98
4372
738
356
344
277
268
224
35
12
0.64
0.63
0.17
-0.12
-0.14
0.07
0.48
-0.06
0.12
hac
om
quarter
10
0.53
0.51
0.13
0.73
0.41
0.52
2.98
4372
738
356
344
277
268
224
20
15
12
0.64
0.63
0.17
-0.12
-0.14
0.07
0.48
0.12
0.10
0.12
hac
om
quarter
15
0.56
0.51
0.09
0.79
0.47
0.59
4.92
4210
738
268
231
224
203
183
173
113
102
74
60
20
15
12
0.63
0.62
0.02
0.01
0.47
-0.12
0.12
0.75
0.17
0.04
-0.06
0.05
0.07
0.01
0.09
hac
lcs
quarter
8
0.51
0.51
0.06
0.86
0.47
0.55
7.31
3327
1807
742
237
202
130
91
90
1.00
-0.17
0.38
0.21
0.37
-0.03
-0.03
-0.03
pam
lcs
quarter
11
0
0.57
0.51
0.05
0.92
0.45
0.56
10.49
4349
742
294
251
201
191
146
126
125
117
84
0.62
0.43
-0.07
-0.16
0.63
0.47
-0.08
0.44
0.67
0.39
0.01
pam
lcs
quarter
15
1
0.57
0.51
0.03
0.94
0.50
0.60
17.86
4055
695
299
224
193
187
172
134
111
107
99
91
89
85
85
0.65
0.50
0.27
-0.26
0.02
0.47
0.15
0.08
0.20
0.56
0.55
1.00
0.27
-0.08
0.00
hac
om
quarter
6
0.50
0.50
0.17
0.68
0.34
0.43
2.00
4996
738
356
277
224
35
0.54
0.65
0.23
-0.10
0.51
-0.06
hac
om
quarter
14
0.56
0.50
0.09
0.79
0.46
0.58
4.92
4210
738
344
268
224
203
183
173
102
74
60
20
15
12
0.63
0.62
-0.15
0.03
0.47
-0.12
0.12
0.75
0.04
-0.06
0.05
0.07
0.01
0.12
hac
lcs
quarter
7
0.49
0.50
0.08
0.81
0.45
0.52
4.96
3327
2044
742
202
130
91
90
1.00
-0.17
0.38
0.37
-0.01
0.07
-0.01
pam
lcs
quarter
4
0
0.51
0.50
0.17
0.79
0.27
0.37
2.37
5250
745
379
252
0.56
0.50
0.00
0.04
pam
lcs
quarter
10
0
0.56
0.50
0.06
0.90
0.43
0.53
8.40
4386
744
318
251
201
190
164
127
125
120
0.61
0.43
-0.10
-0.16
0.64
0.49
-0.07
0.44
0.67
0.35
pam
lcs
quarter
15
0
0.57
0.50
0.03
0.93
0.50
0.61
17.67
4113
736
290
227
201
191
154
146
115
91
85
83
80
79
35
0.63
0.43
-0.06
-0.12
0.61
0.44
0.11
-0.10
0.27
1.00
0.35
0.01
1.00
0.74
0.05
pam
lcs
quarter
4
1
0.53
0.50
0.15
0.82
0.27
0.40
2.90
5340
734
350
202
0.55
0.53
0.11
-0.06
pam
lcs
quarter
3
0
0.45
0.49
0.21
0.74
0.23
0.31
1.59
5502
745
379
0.51
0.52
0.03
pam
lcs
quarter
5
0
0.51
0.49
0.17
0.80
0.31
0.42
2.40
5016
743
349
263
255
0.57
0.49
0.23
0.02
-0.03
pam
lcs
quarter
9
0
0.57
0.49
0.06
0.90
0.42
0.53
8.55
4485
744
318
251
201
190
164
155
118
0.60
0.44
-0.10
-0.10
0.64
0.50
-0.07
0.37
0.37
pam
lcs
quarter
14
0
0.56
0.49
0.04
0.93
0.49
0.58
13.02
4113
742
294
251
201
191
154
146
115
91
85
84
80
79
0.63
0.39
-0.08
-0.23
0.61
0.44
0.11
-0.09
0.27
1.00
0.35
0.01
1.00
0.74
pam
lcs
quarter
3
1
0.45
0.49
0.21
0.74
0.23
0.31
1.59
5502
739
385
0.51
0.53
0.02
pam
lcs
quarter
5
1
0.51
0.49
0.17
0.80
0.31
0.41
2.40
5016
739
349
281
241
0.56
0.50
0.23
-0.05
0.04
hac
lcs
quarter
5
0.48
0.48
0.09
0.78
0.39
0.46
4.16
3529
2044
742
221
90
0.90
-0.13
0.40
-0.04
0.00
hac
lcs
quarter
6
0.48
0.48
0.09
0.79
0.41
0.50
4.21
3529
2044
742
130
91
90
0.90
-0.13
0.40
0.00
0.07
0.00
pam
lcs
quarter
8
0
0.56
0.48
0.08
0.89
0.40
0.51
6.23
4604
744
322
251
201
190
159
155
0.56
0.46
-0.10
-0.07
0.65
0.51
-0.07
0.40
hac
lcs
quarter
4
0.45
0.47
0.12
0.75
0.37
0.42
2.81
3529
2134
742
221
0.90
-0.17
0.40
-0.03
pam
lcs
quarter
7
0
0.56
0.46
0.10
0.87
0.38
0.49
4.87
4758
744
322
240
204
193
165
0.52
0.48
-0.08
0.02
0.65
0.53
-0.07
hac
lcs
quarter
3
0.43
0.44
0.19
0.70
0.25
0.33
1.58
4271
2134
221
0.71
-0.06
-0.02
pam
lcs
quarter
2
0
0.22
0.44
0.44
0.50
0.13
0.14
0.25
5878
748
0.42
0.55
pam
lcs
quarter
2
1
0.22
0.44
0.44
0.50
0.13
0.14
0.25
5878
748
0.42
0.55
hac
lcs
quarter
2
0.36
0.43
0.25
0.63
0.20
0.20
0.91
4271
2355
0.75
-0.14
pam
lcs
quarter
6
0
0.55
0.42
0.13
0.84
0.34
0.46
3.55
4949
744
322
251
201
159
0.47
0.48
-0.08
0.01
0.67
-0.06
Tiempo que demora esta sección: 0 minutos
Para validar la robustez de esta tipología, se implementó un procedimiento de bootstrap paramétrico con 1000 réplicas, comparando la calidad de la solución observada con la obtenida al aplicar el mismo procedimiento de clustering a datos generados bajo un modelo nulo combinado. Este modelo nulo evalúa la calidad de la agrupación en ausencia de estructura real, considerando aspectos combinados de duración y secuencia (comb), así como la secuencia por sí sola (seq), según la metodología propuesta por Studer (2021). Se busca medir cuánta calidad obtenida por la tipología sobrepasa la que habría sido obtenida para datos sin una estructura de conglomerados. El criterio Max T obtiene los máximos puntajes para las trayectorias que asumen estructuras aleatorias sobre las que comparar la solución de conglomerados obtenida.@studer_validating_2021
cat("Hacemos clasificación de pertenencia cluster a las soluciones candidatas y añadimos etiquetas\n")
Hacemos clasificación de pertenencia cluster a las soluciones candidatas y añadimos etiquetas
Código
ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4 <-factor(pamRange_quarter_om$clustering$cluster4,levels=rev(attr( sort(table(pamRange_quarter_om$clustering$cluster4)), "name")),labels=c("6623, Un trimestre, TSM(4)", "6612, Un trimestre, TUS(3)", "6522, Un semestre TSM(1)", "6574, Comorbilidad un trimestre(2)"))ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4_2 <-factor(pamRange_quarter_om2$clustering$cluster4,levels=rev(attr( sort(table(pamRange_quarter_om2$clustering$cluster4)), "name")),labels=c("1, Un trimestre, TSM(1)", "15, Un trimestre, TUS(2)", "32, Un semestre TSM(4)", "49, Comorbilidad un trimestre(3)"))ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om6 <-#pamRange_quarter_om$clustering$cluster6factor(pamRange_quarter_om$clustering$cluster6,levels=rev(attr( sort(table(pamRange_quarter_om$clustering$cluster6)), "name")),labels=c("6623, Un trimestre, TSM(5)", "6612, Un trimestre, TUS(4)", "6522, Un semestre TSM(2)", "6624, TSM, 1 año después, otras causas(6)","6574, Comorbilidad un trimestre(3)", "6268, TSM, 1 año después, TSM(1)" ))ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2_2 <-factor(pamRange_quarter_om2$clustering$cluster2,levels=rev(attr( sort(table(pamRange_quarter_om2$clustering$cluster2)), "name")), labels=c("1, Un trimestre, TSM(1)", "15, Un trimestre, TUS(2)"))ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om3 <-factor(pamRange_quarter_om$clustering$cluster3,levels=rev(attr( sort(table(pamRange_quarter_om$clustering$cluster3)), "name")), labels=c("6623, Un trimestre, TSM(3)", "6612, Un trimestre, TUS(2)", "6522, Un semestre TSM(1)"))ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2 <-factor(pamRange_quarter_om$clustering$cluster2,levels=rev(attr( sort(table(pamRange_quarter_om$clustering$cluster2)), "name")), labels=c("6623, Un trimestre, TSM(2)", "6612, Un trimestre, TUS(1)"))ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_hc_om2 <-factor(om_dist_quarter_c$clustering$cluster2,levels=rev(attr( sort(table(om_dist_quarter_c$clustering$cluster2)), "name")), labels=c("1, Un trimestre, TSM(1)", "2, Un trimestre, TUS(2)"))#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:cat("(==============================================================)\n")
cat("Creamos valores ASW para las soluciones candidatas\n")
Creamos valores ASW para las soluciones candidatas
Código
# Creamos un vector con las columnas llenando con NA si faltan valores# sil_pam_om_clus6_q <-wcSilhouetteObs(as.dist(dist_quarter_om), pamRange_quarter_om$clustering$cluster6, measure="ASW")sil_pam_om_clus4_q <-wcSilhouetteObs(as.dist(dist_quarter_om), pamRange_quarter_om$clustering$cluster4, measure="ASW")sil_pam_om_clus3_q <-wcSilhouetteObs(as.dist(dist_quarter_om), pamRange_quarter_om$clustering$cluster3, measure="ASW")sil_pam_om_clus2_q <-wcSilhouetteObs(as.dist(dist_quarter_om), pamRange_quarter_om$clustering$cluster2, measure="ASW")sil_hc_om_clus2_q <-wcSilhouetteObs(as.dist(dist_quarter_om), om_dist_quarter_c$clustering$cluster4, measure="ASW")#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#
cat("Visualizamos las soluciones\n")seq_plot_pam_om4_q <-ggseqiplot(States_Wide.seq_quarter_t_prim_adm_cens, group= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4,facet_ncol=2, facet_nrow=2, sortv=sil_pam_om_clus4_q) +theme(legend.position ="none")+labs(x="Trimestres", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_pam_om4_q
Visualizamos las soluciones
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.2 minutos
Código
seq_plot_pam_om6_q <-ggseqiplot(States_Wide.seq_quarter_t_prim_adm_cens, group= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om6,facet_ncol=2, facet_nrow=3, sortv=sil_pam_om_clus6_q) +theme(legend.position ="none")+labs(x="Trimestres", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_pam_om6_q
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.2 minutos
Código
seq_plot_pam_om3_q <-ggseqiplot(States_Wide.seq_quarter_t_prim_adm_cens, group= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om3,facet_ncol=2, facet_nrow=2, sortv= sil_pam_om_clus3_q) +theme(legend.position ="none")+labs(x="Trimestres", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_pam_om3_q
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.2 minutos
Código
seq_plot_pam_om2_q <-ggseqiplot(States_Wide.seq_quarter_t_prim_adm_cens, group= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2,facet_ncol=2, facet_nrow=1, sortv= sil_pam_om_clus2_q) +theme(legend.position ="none")+labs(x="Trimestres", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_pam_om2_q
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.2 minutos
La solución que obtuvo mejores ínices de ajuste ASW fue la de 6 conglomerados, mediante el algoritmo PAM y emparejamiento óptimo (OM). No obstante, los conglomerados 6624, TSM, 1 año después, otras causas y 6268, TSM, 1 año después, TSM, obtuvieron valores ASW negativos, por lo que dicha solución se descarta como alternativa. Por tanto, la solución que sigue en términos de valores ASW más altos es la de 4 conglomerados obtenidos mediante el algoritmo PAM y mediante emparejamiento óptimo (OM), con 0,58 vs. en ASW vs. 0,59 de la solución de 6 conglomerados. Las soluciones que le siguen presentan valores de 0,55 en ASW. La solución de 3 conglomerados distingue entre primer ingreso por TSM, TUS y una estadía de un semestre por TSM, pero la trayectoria de estadía por comorbilidad no emerge. Las soluciones de 2 conglomerados obtenidas mediante el algoritmo jerárquico fueron similares que las obtenidas mediante PAM, diferenciando por el motivo de ingreso entre ingresos por TSM y TUS.
Se descartó la solución de 2 conglomerados por ser demasiado simple y, por lo tanto, carecer de valor explicativo en relación con los estados discretos de interés (a saber, las causas de hospitalización). En contraste con una estructura de trayectoria de duración y secuencia aleatoria, los valores obtenidos en los índices Gamma de Hubbert (HG) y Correlación Punto-Biserial (PBC) se encuentran dentro del intervalo esperado para la solución de 7 conglomerados, mientras que para la solución de 4 conglomerados, el índice C de Hubbert también se encuentra dentro de ese intervalo sumado a los otros. Esto significa que sólo el ASW de la solución de 4 conglomerados es mayor que lo esperado. Ahora, en comparación a una estructura de trayectorias aleatoria en términos secuencias, ambas soluciones tienen valores de ASW, HC y PBC superiores al esperable, mientras que en el caso del índice HG sólo es superior al esperable para la solución de 7 conglomerados.
Si observamos la tabla anterior, podemos notar que la solución otbenida mediante el algoritmo jerárquico (HAC) y emparejamiento óptimo (OM) con 2 conglomerados, obtuvo valores de ASW, HG, PBC y HC acorde a los esperados para una estructura de trayectorias aleatorias en términos de duración y secuencia en los índices ASW. Para los índices HG y PBC obtuvo valores más bajos que los esperados, mientras que para el índice HC obtuvo un valor más alto. Si observamos los valores esperados para uan estructura de secuencia aleatoria, la solución de 2 conglomerados obtuvo valores de ASW superiores a los esperados, pero inferiores a los esperados en el índice HG, y dentro de lo esperado en los índices PBC y HC.
En el caso de la solución de 6 conglomerados vs. una estructura aleatoria en términos de duración y secuencia, obtuvo valores ASW superiores a los esperados y más bajos en HC, aunque bajo lo esperado en los índices HG y PBC. Posteriormente, comparándolo con seuencias aleatorias en términos de secuencia, la solución de 6 conglomerados obtuvo valores de ASW, HG, PBC y superiores a los esperados y un índice HC más que para una estructura de trayectorias aleatorias.
Para las soluciones de 2 y 3 conglomerados vs. una estructura de duración y secuencia aleatorias, los valores ASW se encontraban dentro de lo esperado, mientras que los valores de HG, PBC y HC eran inferiores a lo esperado. En el caso de las soluciones de 2 y 3 conglomerados vs. una estructura de secuencia aleatoria, los valores de ASW fueron superiores al esperado e inferiores al esperado en el índice HG para ambas soluciones. No obstante, los valores PBC fueron superiores y los valores HC fueron inferiores para la solución de 3 conglomerados.
Para la solución de 4 conglomerados vs una estructura de secuencia y duración aleatorias, los valores ASW fueron superiores y los HC inferiores al esperado, mientras que los valores HG y PBC fueron inferiores al esperado. Al compararse con una estructura de secuencia aleatoria, los valores ASW y PBC fueron superiores y el valor HC inferiores al esperado, aunque el valor en el índice HG fue más bajo que el esperado.
A continuación se muestra con más detalle el resultado de pruebas de validación mediante bootstraps de la solución obtenida con el algoritmo PAM y emparejamiento óptimo.
Código
# results: hide# fig.show: hideratio_plot=5asw_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_asw, width =800*ratio_plot, height =600*ratio_plot, res=500)hc_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_hc, width =800*ratio_plot, height =600*ratio_plot, res=500)hg_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_hg, width =800*ratio_plot, height =600*ratio_plot, res=500)pbc_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_pbc, width =800*ratio_plot, height =600*ratio_plot, res=500)final_plot_comb <-plot_grid( asw_grob, hc_grob, hg_grob, pbc_grob,ncol =2, # Número de columnasnrow =2, # Número de filasrel_widths =c(1, 1), # Ancho relativo de los gráficosrel_heights =c(1, 1),labels =c("A", "B", "C", "D"), # Etiquetas opcionaleslabel_size =15, # Tamaño de las etiquetasalign ="v", # Alineación vertical de los gráficosaxis ="tb"# Alineación de ejes superior e inferior)ggdraw() +draw_plot(final_plot_comb, x =0, y =0.1, width =1, height =0.9) +draw_text("Área gris: índices de agrupaciones aleatorias; línea negra: índices obtenidos", x =0.05, y =0.05, hjust =0, size =8, lineheight = .8)
Indicadores de calidad vs. bootstrap con secuencias y duraciones aleatorias
cqi_month<-rbind.data.frame(cbind.data.frame(algo="hac", type="om", time="month", k=2:15, corr=F, om_dist_month_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="hac", type="lcs", time="month", k=2:15, corr=F, lcs_dist_month_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="month", k=2:15, corr=F, pamRange_month_om$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="month", k=2:15, corr=T, pamRange_month_om2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="month", k=2:15, corr=F, pamRange_month_lcs$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="month", k=2:15, corr=T, pamRange_month_lcs2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2)))) |> dplyr::select(algo, type, time, k, corr, PBC, ASW, HC, HG, R2, R2sq)tabs_month_clus_sol<-rbind.data.frame(func_tab_range_clus(om_dist_month_c),func_tab_range_clus(lcs_dist_month_c),func_tab_range_clus(pamRange_month_om),func_tab_range_clus(pamRange_month_om2),func_tab_range_clus(pamRange_month_lcs),func_tab_range_clus(pamRange_month_lcs2))#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:# Inicializamos una lista para almacenar los resultadosresultados_list <-list()# Definimos un rango para los clusters a evaluarcluster_range <-2:15# Definimos los métodos y sus variablesmetodos <-list(hac_om =list(data = om_dist_month_c, dist = dist_month_om),hac_lcs =list(data = lcs_dist_month_c, dist = dist_month_lcs),pam_om0 =list(data = pamRange_month_om, dist = dist_month_om),pam_om1 =list(data = pamRange_month_om2, dist = dist_month_om),pam_lcs0 =list(data = pamRange_month_lcs, dist = dist_month_lcs),pam_lcs1 =list(data = pamRange_month_lcs2, dist = dist_month_lcs))# Número máximo de clusters para definir las columnasmax_clusters <-max(cluster_range)# Iteramos sobre cada métodofor (metodo innames(metodos)) {# Creamos un data frame temporal para cada método metodo_result <-data.frame()# Iteramos sobre cada cluster en el rangofor (cluster in cluster_range) {# Construimos el nombre del cluster dinámicamente cluster_name <-paste0("cluster", cluster)# Intentamos calcular los valores de silhouette silhouette_values <-tryCatch(round(summary(silhouette(as.integer(metodos[[metodo]]$data$clustering[[cluster_name]]), as.dist(metodos[[metodo]]$dist)))$clus.avg.widths[attr(rev(sort(table(metodos[[metodo]]$data$clustering[[cluster_name]]))),"names")], 2),error =function(e) rep(NA, cluster) )# Creamos un vector con las columnas llenando con NA si faltan valores silhouette_full <-c(silhouette_values, rep(NA, max_clusters -length(silhouette_values)))# Creamos un data frame temporal con los resultados para este cluster cluster_result <-data.frame(Metodo = metodo,Cluster = cluster,t(silhouette_full) # Transponemos los valores para que cada uno sea una columna )# Nombramos dinámicamente las columnas de silhouettecolnames(cluster_result)[3:(3+ max_clusters -1)] <-paste0("asw", 1:max_clusters)# Añadimos el resultado del cluster al data frame del método metodo_result <-rbind(metodo_result, cluster_result) }# Agregamos los resultados del método a la lista general resultados_list[[metodo]] <- metodo_result}# Combinamos todos los resultados en un único data frameavs_por_cluster_month <-do.call(rbind, resultados_list)# Ordenamos por Método y Clusteravs_por_cluster_month <- avs_por_cluster_month[order(avs_por_cluster_month$Metodo, avs_por_cluster_month$Cluster), ]bind_cols(cqi_month, tabs_month_clus_sol)|> dplyr::mutate(corr= dplyr::case_when(corr==TRUE& algo!="hac"~"1",corr==FALSE& algo!="hac"~"0",T~""), key=paste0(algo,"_",type,corr,"_",k))|>left_join(dplyr::mutate(avs_por_cluster_month, key=paste0(Metodo,"_",Cluster)), by="key")|> dplyr::select(-Metodo, -Cluster) |>`rownames<-`(NULL) |> dplyr::mutate(calc=round(PBC*(1/HC)*HG,2)) |> dplyr::arrange(desc(ASW)) |> dplyr::select(c("algo", "type", "time", "k", "corr", "PBC", "ASW", "HC", "HG", "R2", "R2sq", "calc", paste0("X",1:15), paste0("asw",1:15)))|> (\(df) {assign("asw_month_qci", dplyr::select(df, -"time"), envir = .GlobalEnv) rio::export(df, "_output/sol_conglomerados_tab_month_25.xlsx") knitr::kable(df, "markdown", caption ="CQIs y frecuencias en conglomerados (mensual)") })()
CQIs y frecuencias en conglomerados (mensual)
algo
type
time
k
corr
PBC
ASW
HC
HG
R2
R2sq
calc
X1
X2
X3
X4
X5
X6
X7
X8
X9
X10
X11
X12
X13
X14
X15
asw1
asw2
asw3
asw4
asw5
asw6
asw7
asw8
asw9
asw10
asw11
asw12
asw13
asw14
asw15
hac
om
month
2
0.49
0.86
0.15
0.99
0.03
0.12
3.23
6598
28
0.87
0.06
hac
lcs
month
2
0.58
0.63
0.10
0.91
0.15
0.22
5.28
5740
886
0.76
-0.21
hac
lcs
month
3
0.60
0.61
0.09
0.91
0.18
0.34
6.07
5740
840
46
0.73
-0.16
-0.07
pam
om
month
15
0
0.44
0.45
0.11
0.75
0.39
0.50
3.00
4208
749
585
197
157
117
113
102
93
87
85
67
44
16
6
0.57
0.47
0.33
0.19
0.17
0.12
-0.06
-0.14
-0.07
-0.11
-0.15
-0.18
0.06
0.09
0.14
pam
om
month
14
1
0.44
0.45
0.11
0.74
0.38
0.49
2.96
4208
749
585
197
168
144
116
103
93
88
78
73
18
6
0.57
0.47
0.33
0.19
0.19
-0.03
-0.06
-0.16
-0.07
-0.11
-0.12
-0.19
0.10
0.15
pam
om
month
15
1
0.44
0.45
0.10
0.76
0.39
0.50
3.34
4127
749
585
197
168
144
114
102
93
87
85
84
67
18
6
0.60
0.46
0.31
0.18
0.19
-0.04
-0.06
-0.15
-0.09
-0.12
-0.15
-0.15
-0.19
0.10
0.15
pam
om
month
8
0
0.42
0.44
0.15
0.68
0.33
0.37
1.90
4477
753
620
235
198
140
106
97
0.52
0.48
0.29
0.10
0.22
-0.20
-0.15
-0.25
pam
om
month
9
0
0.41
0.44
0.16
0.68
0.34
0.42
1.74
4450
749
622
198
180
148
135
123
21
0.53
0.49
0.27
0.23
0.02
0.05
-0.20
-0.09
0.02
pam
om
month
10
0
0.42
0.44
0.14
0.70
0.35
0.44
2.10
4408
750
605
197
178
148
119
119
89
13
0.54
0.48
0.30
0.22
0.10
-0.02
-0.06
-0.17
-0.22
0.14
pam
om
month
11
0
0.43
0.44
0.13
0.71
0.36
0.45
2.35
4319
750
605
197
176
148
118
118
95
87
13
0.56
0.47
0.29
0.21
0.13
-0.03
-0.06
-0.17
-0.11
-0.22
0.14
pam
om
month
12
0
0.43
0.44
0.12
0.72
0.37
0.46
2.58
4294
749
585
197
170
148
114
103
93
87
73
13
0.55
0.48
0.34
0.21
0.15
-0.05
-0.05
-0.14
-0.05
-0.18
-0.23
0.14
pam
om
month
13
0
0.43
0.44
0.12
0.72
0.37
0.47
2.58
4294
749
585
197
157
128
115
102
93
86
71
38
11
0.55
0.48
0.34
0.21
0.20
0.11
-0.05
-0.12
-0.05
-0.18
-0.22
-0.06
0.13
pam
om
month
14
0
0.43
0.44
0.12
0.72
0.38
0.49
2.58
4294
749
585
197
157
117
113
103
93
85
67
44
16
6
0.55
0.48
0.34
0.21
0.17
0.12
-0.05
-0.14
-0.05
-0.14
-0.18
0.06
0.09
0.14
pam
om
month
7
1
0.39
0.44
0.17
0.64
0.31
0.40
1.47
4539
754
729
239
198
145
22
0.52
0.49
0.24
-0.08
0.24
-0.20
0.03
pam
om
month
8
1
0.42
0.44
0.15
0.68
0.33
0.37
1.90
4477
753
620
235
198
140
106
97
0.52
0.48
0.29
0.10
0.22
-0.20
-0.15
-0.25
pam
om
month
9
1
0.42
0.44
0.15
0.68
0.34
0.43
1.90
4477
753
620
235
198
139
99
93
12
0.52
0.48
0.29
0.10
0.22
-0.19
-0.02
-0.22
0.15
pam
om
month
10
1
0.42
0.44
0.14
0.70
0.35
0.44
2.10
4408
750
605
197
181
148
119
116
89
13
0.54
0.48
0.30
0.22
0.07
-0.01
-0.17
-0.04
-0.22
0.14
pam
om
month
11
1
0.43
0.44
0.13
0.71
0.36
0.45
2.35
4319
750
605
197
179
148
118
115
95
87
13
0.56
0.47
0.29
0.21
0.10
-0.02
-0.17
-0.04
-0.11
-0.22
0.14
pam
om
month
12
1
0.43
0.44
0.12
0.72
0.37
0.46
2.58
4294
749
585
197
169
148
112
111
93
78
77
13
0.55
0.48
0.34
0.21
0.19
-0.07
-0.15
-0.05
-0.05
-0.22
-0.18
0.14
pam
om
month
13
1
0.43
0.44
0.12
0.72
0.37
0.49
2.58
4294
749
585
197
169
144
112
111
93
75
73
18
6
0.55
0.48
0.34
0.21
0.19
-0.03
-0.15
-0.05
-0.05
-0.12
-0.19
0.10
0.15
hac
om
month
10
0.39
0.43
0.17
0.64
0.29
0.42
1.47
4312
995
743
314
186
28
25
18
3
2
0.60
-0.17
0.51
-0.16
0.35
-0.17
-0.02
-0.07
0.49
0.75
hac
om
month
11
0.40
0.43
0.15
0.66
0.31
0.45
1.76
4312
743
728
314
267
186
28
25
18
3
2
0.57
0.51
0.08
-0.19
-0.28
0.34
-0.18
-0.03
-0.09
0.49
0.75
hac
om
month
12
0.40
0.43
0.15
0.66
0.32
0.46
1.76
4312
743
728
314
267
186
28
18
15
10
3
2
0.57
0.51
0.08
-0.19
-0.28
0.34
-0.18
-0.09
-0.04
0.39
0.44
0.75
pam
om
month
5
0
0.35
0.43
0.22
0.59
0.27
0.29
0.94
4645
760
759
256
206
0.50
0.20
0.46
-0.18
0.21
pam
om
month
6
0
0.39
0.43
0.18
0.63
0.30
0.35
1.36
4590
754
724
229
205
124
0.50
0.48
0.26
-0.06
0.21
-0.25
pam
om
month
7
0
0.41
0.43
0.16
0.67
0.32
0.36
1.72
4492
751
705
218
203
140
117
0.52
0.47
0.29
-0.06
0.21
-0.12
-0.26
pam
om
month
6
1
0.39
0.43
0.18
0.63
0.30
0.34
1.36
4590
754
633
320
205
124
0.50
0.48
0.26
-0.03
0.21
-0.25
pam
lcs
month
2
0
0.25
0.43
0.39
0.48
0.10
0.12
0.31
5646
980
0.52
-0.05
pam
lcs
month
2
1
0.25
0.43
0.39
0.48
0.10
0.12
0.31
5646
980
0.52
-0.05
hac
om
month
9
0.36
0.42
0.19
0.62
0.28
0.40
1.17
4340
995
743
314
186
25
18
3
2
0.58
-0.16
0.52
-0.15
0.36
-0.02
-0.07
0.49
0.75
hac
om
month
13
0.40
0.42
0.15
0.67
0.33
0.48
1.79
4312
743
691
314
267
186
37
28
18
15
10
3
2
0.53
0.51
0.22
-0.20
-0.29
0.34
-0.11
-0.18
-0.09
-0.08
0.39
0.44
0.75
hac
om
month
14
0.40
0.42
0.15
0.67
0.34
0.49
1.79
4312
743
691
269
267
186
45
37
28
18
15
10
3
2
0.53
0.51
0.22
-0.18
-0.29
0.34
0.24
-0.12
-0.18
-0.10
-0.08
0.32
0.44
0.75
hac
om
month
15
0.40
0.42
0.15
0.67
0.34
0.50
1.79
4312
743
691
269
267
186
45
37
28
15
15
10
3
3
2
0.53
0.51
0.22
-0.18
-0.29
0.34
0.24
-0.12
-0.18
-0.04
-0.08
0.32
0.43
0.19
0.75
hac
lcs
month
15
0.37
0.42
0.05
0.84
0.45
0.59
6.22
2786
1412
726
641
435
176
175
145
60
24
16
13
12
3
2
1.00
-0.22
0.28
0.38
-0.23
-0.16
0.12
0.00
-0.19
-0.02
-0.01
0.07
0.33
0.46
0.75
pam
om
month
3
0
0.30
0.42
0.28
0.51
0.21
0.21
0.55
4874
987
765
0.50
-0.06
0.47
pam
om
month
3
1
0.30
0.42
0.28
0.51
0.21
0.21
0.55
4874
987
765
0.50
-0.06
0.47
pam
om
month
5
1
0.35
0.42
0.22
0.59
0.27
0.29
0.94
4645
759
666
350
206
0.50
0.46
0.19
-0.15
0.21
pam
lcs
month
10
1
0.42
0.42
0.11
0.78
0.36
0.49
2.98
4269
737
696
252
173
133
128
112
105
21
0.54
0.37
0.30
-0.03
0.20
-0.20
-0.02
0.02
-0.16
0.03
pam
lcs
month
14
1
0.44
0.42
0.09
0.83
0.40
0.53
4.06
4177
714
573
174
164
148
128
124
104
97
88
60
53
22
0.55
0.43
0.29
0.17
0.48
-0.14
-0.04
-0.18
-0.15
0.04
-0.15
-0.14
-0.14
0.01
pam
lcs
month
15
1
0.44
0.42
0.08
0.84
0.41
0.54
4.62
4143
701
547
168
162
150
146
119
107
85
84
80
60
52
22
0.55
0.46
0.34
0.20
0.48
-0.26
-0.13
-0.01
-0.11
0.10
-0.15
-0.09
-0.14
-0.16
0.00
hac
om
month
3
0.25
0.41
0.40
0.42
0.14
0.22
0.26
5667
931
28
0.41
0.39
0.04
hac
lcs
month
4
0.49
0.41
0.12
0.80
0.26
0.38
3.27
5014
840
726
46
0.50
-0.18
0.47
-0.07
hac
lcs
month
5
0.49
0.41
0.12
0.80
0.26
0.41
3.27
5014
840
726
31
15
0.50
-0.18
0.47
-0.07
0.21
hac
lcs
month
6
0.49
0.41
0.12
0.80
0.27
0.44
3.27
5014
816
726
31
24
15
0.50
-0.17
0.47
-0.08
0.02
0.21
pam
lcs
month
8
0
0.40
0.41
0.14
0.73
0.33
0.46
2.09
4395
749
727
254
184
150
142
25
0.52
0.36
0.22
0.01
0.13
-0.12
-0.18
0.02
pam
lcs
month
9
0
0.42
0.41
0.13
0.76
0.35
0.48
2.46
4326
749
700
228
182
169
132
116
24
0.53
0.34
0.28
0.04
0.14
-0.10
-0.18
-0.15
0.03
pam
lcs
month
10
0
0.42
0.41
0.12
0.78
0.36
0.49
2.73
4261
746
693
197
180
171
127
119
108
24
0.54
0.34
0.28
0.06
0.15
-0.12
-0.19
-0.16
0.03
0.03
pam
lcs
month
11
0
0.42
0.41
0.12
0.78
0.37
0.50
2.73
4278
746
599
180
170
164
128
114
112
111
24
0.53
0.34
0.19
0.15
-0.12
0.51
-0.19
-0.15
0.02
-0.08
0.01
pam
lcs
month
13
0
0.44
0.41
0.09
0.81
0.39
0.52
3.96
4217
735
578
178
164
164
112
111
101
101
91
49
25
0.54
0.37
0.28
0.14
0.48
-0.10
-0.09
-0.15
0.03
-0.15
-0.15
-0.12
0.00
pam
lcs
month
15
0
0.44
0.41
0.08
0.83
0.41
0.54
4.56
4161
727
563
190
181
116
109
106
104
91
87
62
60
49
20
0.54
0.40
0.26
-0.11
0.13
-0.16
-0.01
-0.16
0.79
-0.15
0.18
-0.14
-0.06
-0.13
0.00
pam
lcs
month
8
1
0.42
0.41
0.13
0.76
0.34
0.43
2.46
4342
742
709
261
176
145
137
114
0.53
0.37
0.28
-0.02
0.20
-0.19
-0.03
-0.18
pam
lcs
month
11
1
0.42
0.41
0.12
0.78
0.37
0.50
2.73
4285
741
602
172
168
152
140
134
116
101
15
0.53
0.35
0.22
0.21
0.48
-0.13
-0.07
-0.20
0.01
-0.15
0.10
pam
lcs
month
12
1
0.43
0.41
0.11
0.79
0.38
0.51
3.09
4267
739
580
172
164
148
134
128
103
98
71
22
0.52
0.36
0.27
0.19
0.49
-0.14
-0.04
-0.20
0.02
-0.10
-0.18
0.02
pam
lcs
month
13
1
0.44
0.41
0.09
0.81
0.39
0.52
3.96
4217
732
578
177
167
148
134
109
99
96
92
55
22
0.54
0.37
0.28
0.16
0.47
-0.14
-0.05
-0.15
-0.15
0.04
-0.15
-0.15
0.02
hac
om
month
8
0.36
0.40
0.19
0.61
0.26
0.38
1.16
4340
995
929
314
25
18
3
2
0.59
-0.15
0.31
-0.14
-0.02
-0.06
0.49
0.77
pam
om
month
2
0
0.16
0.40
0.49
0.38
0.11
0.10
0.12
5857
769
0.38
0.52
pam
om
month
4
0
0.32
0.40
0.26
0.55
0.24
0.26
0.68
4837
768
765
256
0.45
0.21
0.46
-0.17
pam
om
month
2
1
0.16
0.40
0.49
0.38
0.11
0.10
0.12
5857
769
0.38
0.52
pam
om
month
4
1
0.32
0.40
0.26
0.54
0.24
0.25
0.66
4837
765
674
350
0.45
0.46
0.21
-0.15
pam
lcs
month
12
0
0.43
0.40
0.11
0.79
0.38
0.51
3.09
4267
744
580
180
165
164
127
112
112
102
48
25
0.52
0.34
0.27
0.14
-0.10
0.49
-0.19
-0.15
-0.09
0.02
-0.13
0.00
pam
lcs
month
14
0
0.44
0.40
0.09
0.82
0.40
0.53
4.01
4201
735
568
178
171
160
110
104
100
98
91
50
38
22
0.53
0.37
0.25
0.13
-0.08
-0.11
-0.15
0.79
-0.15
0.03
-0.15
-0.12
0.09
-0.03
pam
lcs
month
7
1
0.40
0.40
0.15
0.73
0.31
0.39
1.95
4395
744
727
258
178
176
148
0.52
0.37
0.22
-0.01
0.18
-0.22
-0.19
pam
lcs
month
9
1
0.42
0.40
0.13
0.76
0.35
0.43
2.46
4351
745
601
204
176
158
145
133
113
0.51
0.36
0.22
-0.21
0.20
0.53
-0.20
-0.01
-0.19
hac
om
month
4
0.33
0.39
0.22
0.58
0.22
0.27
0.87
4340
1327
931
28
0.60
-0.24
0.29
0.00
hac
om
month
5
0.34
0.39
0.21
0.58
0.23
0.30
0.94
4340
1309
931
28
18
0.60
-0.21
0.29
-0.01
-0.04
hac
om
month
6
0.34
0.39
0.21
0.58
0.23
0.33
0.94
4340
1309
929
28
18
2
0.59
-0.21
0.31
-0.01
-0.04
0.77
hac
om
month
7
0.34
0.39
0.21
0.58
0.24
0.36
0.94
4340
1309
929
25
18
3
2
0.59
-0.21
0.31
0.01
-0.05
0.49
0.77
hac
lcs
month
14
0.37
0.39
0.07
0.79
0.43
0.58
4.18
2786
2053
726
435
176
175
145
60
24
16
13
12
3
2
1.00
-0.17
0.28
-0.20
-0.13
0.12
0.12
-0.16
-0.01
-0.01
0.07
0.33
0.46
0.75
pam
lcs
month
6
0
0.37
0.39
0.18
0.67
0.29
0.36
1.38
4512
753
745
268
185
163
0.48
0.37
0.22
-0.05
0.18
-0.20
pam
lcs
month
7
0
0.38
0.39
0.16
0.69
0.32
0.41
1.64
4467
753
711
258
181
145
111
0.47
0.36
0.30
-0.01
0.20
-0.03
-0.18
hac
lcs
month
11
0.37
0.38
0.08
0.76
0.39
0.55
3.52
2961
2053
726
435
236
145
29
24
12
3
2
0.88
-0.12
0.31
-0.19
-0.19
0.12
-0.08
-0.01
0.37
0.47
0.76
hac
lcs
month
12
0.37
0.38
0.08
0.76
0.40
0.56
3.52
2961
2053
726
435
236
145
24
16
13
12
3
2
0.88
-0.12
0.31
-0.19
-0.19
0.12
-0.01
0.00
0.08
0.33
0.46
0.76
hac
lcs
month
13
0.37
0.38
0.08
0.76
0.40
0.57
3.52
2961
2053
726
435
176
145
60
24
16
13
12
3
2
0.88
-0.12
0.31
-0.19
-0.13
0.12
-0.16
-0.01
-0.01
0.07
0.33
0.46
0.76
hac
lcs
month
7
0.34
0.37
0.11
0.71
0.36
0.46
2.19
2961
2053
816
726
31
24
15
0.88
-0.12
-0.20
0.31
-0.08
0.00
0.21
hac
lcs
month
8
0.34
0.37
0.11
0.71
0.37
0.49
2.19
2961
2053
816
726
31
24
12
3
0.88
-0.12
-0.20
0.31
-0.10
0.00
0.41
0.48
hac
lcs
month
9
0.34
0.37
0.11
0.71
0.37
0.51
2.19
2961
2053
816
726
29
24
12
3
2
0.88
-0.12
-0.20
0.31
-0.05
0.00
0.40
0.47
0.76
hac
lcs
month
10
0.36
0.37
0.10
0.74
0.38
0.53
2.66
2961
2053
726
671
145
29
24
12
3
2
0.88
-0.12
0.31
-0.23
0.13
-0.06
-0.01
0.37
0.47
0.76
pam
lcs
month
5
0
0.36
0.37
0.20
0.64
0.27
0.34
1.15
4697
753
745
268
163
0.43
0.37
0.25
-0.05
-0.19
pam
lcs
month
4
0
0.32
0.36
0.24
0.58
0.23
0.27
0.77
4790
778
760
298
0.42
0.21
0.35
-0.19
pam
lcs
month
3
0
0.26
0.34
0.31
0.50
0.18
0.18
0.42
4885
980
761
0.42
-0.09
0.37
pam
lcs
month
3
1
0.26
0.34
0.31
0.50
0.18
0.18
0.42
4885
980
761
0.42
-0.09
0.37
pam
lcs
month
4
1
0.34
0.30
0.23
0.59
0.25
0.30
0.87
4803
754
676
393
0.31
0.39
0.47
-0.24
pam
lcs
month
5
1
0.36
0.29
0.20
0.64
0.28
0.35
1.15
4697
750
624
392
163
0.30
0.38
0.55
-0.16
-0.20
pam
lcs
month
6
1
0.37
0.25
0.18
0.66
0.30
0.36
1.36
4532
748
536
453
187
170
0.25
0.38
0.71
-0.19
0.14
-0.21
Tiempo que demora esta sección: 0 minutos
Código
frobenius_norm(as.matrix(func_tab_range_clus(pamRange_quarter_lcs)), as.matrix(func_tab_range_clus(pamRange_quarter_lcs2)))invisible("Frobenius norm: ||A-B||_F=\sqrt{\sum{i,j}(A_{ij}-B_ij)^2}")invisible("Por ahora lo dejamos pasar, es para comparar matrices y sus diferencias")
Tiempo que demora esta sección: 0 minutos
La mayoría de las soluciones presentaron conglomerados con tamaños muy pequeños (n < 30), lo que limita su generalización y sugiere una estabilidad presumiblemente baja, o bien con valores ASW negativos, lo que añade evidencia a problemas de estabilidad de la solución.
Código
asw_month_qci$label<-with(asw_month_qci, paste0(algo, "_",type, "_", k, "_", ifelse(corr==1,1,ifelse(corr==0,0,0))))plot(sort(asw_month_qci$ASW, decreasing =TRUE), type ="b", pch =19, xlab ="", xaxt ="n",ylab ="Valor de ASW", main ="", col ="blue")axis(1, at =1:length(asw_month_qci$label), labels = asw_month_qci$label, las =2, cex.axis =0.60)
Gráfico de codo, ASW con etiquetas soluciones de conglomerados (soluciones de res. mensual)
Tiempo que demora esta sección: 0 minutos
Sólo las soluciones de 2 y 3 conglomerados obtuvieron valores ASW en cada conglomerados superiores a 0, es decir, no obtuvieron conglomerados negativos.
cat("Creamos valores ASW para las soluciones candidatas\n")
Creamos valores ASW para las soluciones candidatas
Código
# Creamos un vector con las columnas llenando con NA si faltan valores# sil_hc_om_clus2_m <-wcSilhouetteObs(as.dist(dist_month_om), om_dist_month_c$clustering$cluster2, measure="ASW")sil_hc_om_clus3_m <-wcSilhouetteObs(as.dist(dist_month_om), om_dist_month_c$clustering$cluster3, measure="ASW")sil_pam_om_clus2_m <-wcSilhouetteObs(as.dist(dist_month_om), pamRange_month_om$clustering$cluster2, measure="ASW")#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:
cat("(==============================================================)\n")cat("Visualizamos las soluciones\n")seq_plot_hc_om2_m <-ggseqiplot(States_Wide.seq_month_t_prim_adm_cens, group= ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_hc_om2,#om_dist_month_c$clustering$cluster2,facet_ncol=2, facet_nrow=1, sortv=sil_hc_om_clus2_m) +theme(legend.position ="none")+labs(x="Meses", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_hc_om2_m
(==============================================================)
Visualizamos las soluciones
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.4 minutos
Código
seq_plot_hc_om3_m <-ggseqiplot(States_Wide.seq_month_t_prim_adm_cens, group= ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_hc_om3,facet_ncol=2, facet_nrow=2, sortv=sil_hc_om_clus3_m) +theme(legend.position ="none")+labs(x="Meses", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_hc_om3_m
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.4 minutos
Código
seq_plot_pam_om2_m <-ggseqiplot(States_Wide.seq_month_t_prim_adm_cens, group= ing_dt_ing_calendar_month_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om2,facet_ncol=2, facet_nrow=2, sortv=sil_pam_om_clus2_m) +theme(legend.position ="none")+labs(x="Meses", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_pam_om2_m
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.4 minutos
De las solciones de conglomerados que no obtuvieron algún conglomerado con valores ASW negativos, la solución que obtuvo mejores ínices de ajuste ASW fue la de 2 conglomerados, mediante el algoritmo jerárquico y emparejamiento óptimo (OM) con valores ASW 0,86. Sin embargo, la agrupación es muy poco informativa respecto a los motivos (ej., un episodio de cualquier tipo, vs. múltiples episodios). Por otra parte, la agrupación con múltiples episodios está compuesta por menos de 30 observaciones.
La siguiente solución obtuvo índices de ajuste subóptimos (ASW= 0,41), aunque logró distinguir por diagnóstico del episodio hospitalario, aunque persiste un tercer grupo con múltiples episodios con menos de 30 observaciones.
La tercera y cuarta soluciones obtuvieron un índice de ajuste supótimos (ASW=0,40), aunque distinguir entre ingresos con diagnósticos TSM el primer semestre, vs. el primer trimestre por TSM.
A continuación, vemos los índices de calidad de la solución mediante remuestreos bootstrap con al menos 1,000 replicaciones.
Sólo los índices ASW mostraron mejores índices de calidad que los esperados para una estructura aleatoria en términos de secuencias, mientras que el resto de los índices se encuentran dentro o por debajo de los rangos esperados. Para el resto de los indicadores, se encuentran por debajo de una estructura de secuencias aleatorias en duración y secuencia.
A continuación se muestra el resultado de pruebas de validación mediante bootstraps.
Código
# results: hide# fig.show: hideratio_plot=5asw_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_asw, width =800*ratio_plot, height =600*ratio_plot, res=500)hc_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_hc, width =800*ratio_plot, height =600*ratio_plot, res=500)hg_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_hg, width =800*ratio_plot, height =600*ratio_plot, res=500)pbc_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_pbc, width =800*ratio_plot, height =600*ratio_plot, res=500)final_plot_comb_m <-plot_grid( asw_grob_m, hc_grob_m, hg_grob_m, pbc_grob_m,ncol =2, # Número de columnasnrow =2, # Número de filasrel_widths =c(1, 1), # Ancho relativo de los gráficosrel_heights =c(1, 1),labels =c("A", "B", "C", "D"), # Etiquetas opcionaleslabel_size =15, # Tamaño de las etiquetasalign ="v", # Alineación vertical de los gráficosaxis ="tb"# Alineación de ejes superior e inferior)ggdraw() +draw_plot(final_plot_comb_m, x =0, y =0.1, width =1, height =0.9) +draw_text("Área gris: índices de agrupaciones aleatorias; línea negra: índices obtenidos", x =0.05, y =0.05, hjust =0, size =8, lineheight = .8)
Indicadores de calidad vs. bootstrap con secuencias y duraciones aleatorias
